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Abstract. A detailed analysis of experimentally obtained curvilinear crack path trajectories 
formed in a heterogeneous stress field is presented. The experimental crack path trajectories, 
described in [1], were used as data for the numerical simulations, recreating the actual stress 
field governing the development of the crack path. Thus, the current theories of crack curving 
and kinking could be examined by comparing them with the actual stress field parameters as 
they develop along the experimentally observed crack path. The experimental curvilinear crack 
path trajectories were formed in the tensile specimens with a hole positioned in the vicinity of a 
potential crack path [1]. The numerical simulation, based on the solution of equivalent 
boundary value problems with the possible perturbations of the crack path, is presented here. 


1. Introduction 

A commonly observed curvilinear crack propagation pattern is attributed to a mixed mode 
stress state in the vicinity of the crack tip. Generally, a curvilinear, wave-like, crack path is 
beneficial for material toughening [2], The understanding of the mechanics of crack path 
kinking and curving is important for future design of material compositions capable of 
deflecting the crack propagation trajectory. This work aims to evaluate the crack path 
deflection criteria described in the literature and to state typical aspects of the 
experimentally observed crack path trajectory formation features. The approach 
undertaken here is essentially different from existing work of similar aims. This 
investigation is based strictly on experimental observations which are interpreted through a 
numerical analysis of the observed crack patterns with finite length of the curving segments. 
This makes it different from crack path analysis or predictions based on an infinitesimally 
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small segment of crack kinking and curving, as compared with the length scale of the 
problem [3-16]. 

The experimental data were obtained and reported by Chudnovsky et al. [1]. Briefly, 
the experimental procedure can be described as follows: 

A polystyrene single notch specimen with a hole positioned near the expected crack 
path (Fig. 1 a) was subjected to a low cycle fatigue loading. Cycle loading was used to ensure 
stable crack growth. This particular material does not exhibit a time dependent fracture 
mechanism which could be attributed to fatigue crack propagation; therefore, in the 
described experiment, cycle loading was used to prevent dynamic unstable crack growth. 
Several geometries were used to obtain different intensities of crack hole interaction and to 
ensure that the region of this interaction would be confined to the center of the specimen, 
thus preventing influence from the ends of the specimens. The crack propagation described 
in [1] is accompanied by a damage zone typical for this material, the mechanical nature of 
which is formation of crazes. The typical size of the damage zone is significantly smaller 
than the geometrical parameter of the problem (diameter of the hole). Thus, tw'o 
observations are possible; observation on microscale and on macroscale. The first, which 
was a primary subject of [1], is aimed at processes typical to the region in the vicinity of the 
propagating crack tip, so the size of the damage zone becomes the length parameter of the 
observed physical process. The latter considers the crack with the adjacent propagating 
damage zone as one crack interacting with a hole. This would be the view of an observer 
from a certain distance, when the damage zone appears just as some fuzziness along the 
crack surface. Naturally, the results of the analysis on macroscale can serve as boundary 
value data for the analysis on microscale. 

A summary of the experimental data, observed in [1], is given in Figure 2. Six 
experimental curves, normalized by the hole radius, in the effective region are given. Three 
different cases of the crack interaction with a hole are represented here: crack attraction, 
cases one, six and five; attraction and repulsion, cases three and two; and, finally, no effect 
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on the crack path, case four (the numbering sequence here has no significance; it is left in 
the given form for data reference only). 

2. Analytical Formulation and Numerical Procedure. 

The dimension of the area of crack - hole interaction is significantly smaller than the crack 
length, and smaller than the width of the specimen, as well. The size of this area can be 
judged by the distance from the hole to the first crack path deviation from the rectilinear 
trajectory. Therefore, the effective region can be assumed to be located in the infinite plate 
and to be governed by the applied stress intensity factor. Thus, the small scale model may 
be employed. The applied (or remote) stress intensity factor in this case is the stress 
intensity factor which would act on a straight crack in an equal plate without a hole. The 
length of this crack is equal to a projection of the curvilinear crack onto a horizontal 
direction, Figure 1, a) and b). Thus, the solutions for the local (actual) stress intensity 
factors and other crack growth criteria are obtained and given as in nondimensional form, 
in most cases as a value normalized by the applied stress intensity factor Kf 1 . 

The stress field of a curvilinear crack interacting with a hole is obtained by employing 
the numerical approach developed for a small scale model in [2,18]. The stress field is 
represented by analytic potentials <p(z) and i>(z), Muskhelishvili [17]. The stress tensor 
components in the coordinate system tangential (r ),£ ) to the crack surface can be written as 

° rjr} + ia £r) “ <*>' + <£' + e 2ld {z<p~ (z) + V'(z)), (2.1) 

- i 0 £ V = 4>' + e 21 ^ (z<t>~ (z) + V> ' (z) ) • ( 2 . 2 ) 

where 6 is the angle between the positive direction of x axis and tangential to the crack 
surface £ . 

The crack is represented as an array of dislocations along its curvilinear trajectory. 
Resulting stress potentials are 
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<#>'( z ) = 271 i~t + ^s (b(s ) ,z,t ) ds ' t=t x (s)+it y (s) 

J L L -I (2.3) 

^' (Z) = ~ 271 } L _ f^t “ TF^-^t) 2 " ^ s ( b (s) ,z,t)]ds (2.4) 
Here, the integrals are taken along the crack surface L with 5 , a line length of the crack path 
starting at the crack tip, being chosen as an integration variable. b(s) is an unknown 
dislocation density distribution along the crack, and t is an integration point on the crack 
corresponding to the integration variable s. a here is a standard coefficient, a = E/4(l- v 2 ) 
in plain strain case, and a = E/4 in plain stress case. The singular terms under the integrals 
(2.3) and (2.4) represent singular parts of stress potentials corresponding to a dislocation 
with Burgers vector b(s)ds positioned at t. The analytic functions <p$ and V» s represent 
regular parts of these potentials resulting from the interaction of the dislocation with a hole, 
thus satisfying the traction free boundary conditions on the hole surface. These potentials 
are, [18], 

4 > s (b,z,t,) = Jg(b,z-a,t-a) (2.5) 

A _ A 

^g(b/Z/t) = ^ s (b, z-a, t-a) - a<£" (b, z-a, t-a) (2.6) 

where a is center of the hole, R is the radius of the hole, and potentials corresponding to the 
case a = 0 are 
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density function b(s ). The integral equation is obtained by substitutions of (2.3) and (2.4), 
with the use of (2.5-8), into (2.1) and then satisfying (2.9). The resulting integral equation is 
too long to reproduce in explicit form, and not all the details are important for development 
of the numerical procedure. The resulting equation can be written in the form 


fr N(b(s) ,b(s) ) 

I 0 L t( s ) " z ( v ) 


+ P(b(s) , b (s) , K(s,v) ) 


ds 


0 . 


( 2 . 10 ) 


Functions N and P in (2.10) are linear functions of their arguments. K(s,v) here is a 
Fredholm type kernel of the integral equation (2.10). In reality, functions b(s ) and the 
conjugate to it have different Fredholm kernels not shown in (2.10) explicitly in order to 
emphasize just the essential features of the equation. Important components of it are 
functions t(s ) and z(v). These functions represent the transformation of the integral 
equation in a complex plane into a line integral equation in terms of real variables s and v 
along the crack path L. An important restriction on possible crack trajectories is a 
requirement that 


|t(s) - z(v)| — ♦ const|s - v| as (s - v) — » 0. (2.11) 

Variables ^ and v are the real variables on the integration curve, measured as the curve 
length starting from the crack tip. With condition (2.8) equation (2.7) becomes a first kind 
Cauchy type singular integral equation and the integral is understood in terms of Cauchy 
principal value. 

Following [2], the dislocation density function b(s) in the form defined in (2.3-4), and 
consequently in the equation (2.10), does not correspond to a standard dislocation density 
used in fracture mechanics. The difference is in the coefficient. To use the standard 
definition, ds has to be replaced by dt in equations (2.3,4) and (2.10). This definition was 
used for the convenience of the numerical scheme only, and the difference is accounted for 
in the determination of the stress intensity factors. Essentially, the coefficient dt/ds has been 
absorbed into the unknown dislocation density for computational convenience. However, 
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the restriction on the integration path, following from the requirements of continuity of the 
derivative dt/ds, remains. Thus, the crack trajectory has to be a smooth curve with a 
continuous derivative. 

Equation (2.10) is a homogeneous singular integral equation of the first kind on a 
semi-infinite interval; thus, the stabilization procedure has to be used in the formulation of 
the numerical procedure [18]. The regular numerical scheme is unstable because the 
solution of singular integral equations of this kind, generally, is not unique. The families of 
solutions of this equation are determined by the type of additional conditions imposed on 
the solution of the integral equation. One of these conditions can be an asymptotic behavior 
of the unknown function as 5 - «. The behavior of b(s) for large s should correspond to 
remote loading conditions, or, as in the considered case, equivalently, to the reference 
loading. Thus, the dislocation density function can be written in the following form 


b(s) 


£ 1*1 + i (1 ~ e ~ s2 ) 1 

r 

J s aj 2 7t s -I 


( 2 . 12 ) 


where K m is remote (applied) stress intensity factor and function fi(s) is an unknown 
function bounded on the integration interval and 


/?(s)— - o(s'^) ass — * w . 


(2.13) 


The form chosen here is convenient for computations in the vicinity of the crack tip and the 
evaluation of the stress intensity factors at the crack tip, that is 


Ki + iK 


II 


K 


= i /9(0)e 




(2.14) 


6 0 here is the angle of the tangential at the crack tip. K with superscript o corresponds to 
local values of the stress intensity factors for Mode I and Mode II accordingly. 

After substitution of (2.12) into (2.10) the equation is mapped onto a finite interval and 
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then the collocation procedure is applied, in the form given by Rubinstein, [18], for the 
semi-infinite interval. This numerical procedure is based on the technique introduced by 
Erdogan and Gupta, [19]. 

The collocation scheme uses Gaus - Chebyshev quadrature formula and requires the 
node distribution along the roots of Chebyshev polynomials of the first and second kind. 
To secure the accuracy of the integration procedure, a nonlinear equation 


x k 1 

*(l + (y- (x)) 2 )*dx 

x 


(2 . 15) 


was solved to establish the relation between node and coordinates on the trajectory x,y. 
(Equation (2.13) corresponds to a mapped state.) s £ here is a value of corresponding 
Chebyshev root, ,r° is a coordinate of the crack tip, and function y(x) specifies the crack 
path. 

3. Observation of the experimental trajectories. 

The experimental crack path trajectories were analyzed by employing the technique 
described above. The experimental trajectories were approximated by polynomials using 
the list square method. These trajectories were used as initial data for the computational 
scheme. The computations were conducted for variable crack length along the given 
trajectories and for cases with a possibility of crack path deviations from the given 
experimental path. 

The aim of the analysis is to establish a dominant characteristic parameter determining 
the crack path formation. In pursuit of this, the main fracture parameters are analyzed; 
namely, the focus is on the variations of the stress intensity factors and on the variation of 
the energy release rate due to the crack advance. The experimental specimen is illustrated 
in Figure la); the geometry of the active region is exaggerated there in order to make it 
visible on the figure. The maximal linear dimension of the region, where any significant 
interaction between the hole and a crack take place, is much smaller than the width of the 
specimen; therefore, small scale analysis can be applied. The resulting data are given as 
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nondimensional values which were obtained by normalizing the stress intensity factors by 
the reference value of Aj (illustrated in Figure lb); the energy release rate is given as a 
value normalized by the corresponding reference value. 

3.1 Variations of the crack tip parameters along the trajectories. The results of the 
computations for crack growth along the experimentally obtained trajectories are given in 
Figures 3-5. Cases marked on these figures correspond to trajectories marked on Figure 2. 
The data depicted in Figures 3-5 demonstrate that the crack growth is not characterized by 
monotonic increase of the stress intensity factor Aj or by zero value of Ajj as is often 
expected. Trajectories 3 and 4 are at the distance of approximately 2 R ( R is the radius of 
the hole), and cracks growing along these trajectories experience very weak interaction with 
the hole, although Ajj rises to about 5% of Aj. Remarkably, cracks growing along the 
trajectories 1, 5, and 6 experience a very similar history for all analyzed parameters. 
Trajectory 2 is located in the region R < y < 2R (y is a vertical coordinate of the trajectory 
origin with the center of the hole at y = 0). This region is evidently characterized by a crack 
- hole interaction sufficient to turn the crack towards the hole, but not sufficient to maintain 
this direction, and eventually the crack returns to its original course. 

The history of Ajp Figure 4, is given in values normalized by the reference value of K\, 
which may be misleading in the way that the resulting data become relatively small. The 
values of Ajj vary in the range of ± 15% of the reference Mode I stress intensity factor. The 
generally accepted assumption that cracks grow in the direction of zero Ajj is not supported 
by the obtained data. The possibility of numerical error was ruled out after multiple 
repetitions of the computations with different values of the parameters critical for the 
numerical scheme, such as number of nodes in the integration procedure, degree of the 
approximation polynomials, and so on. The frequency of the registered data points on the 
trajectories is sufficiently higher than the frequency of the oscillation of the values of Ajj, 
which rules out the possibility of insufficient mesh size in trajectory digitizing. Thus, it was 
made sure that the accuracy of the method of analysis is significantly higher than the 
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observed amplitude of K\\ oscillation. It appears that any crack path trajectory with finite 
curvature will have some finite value of Mode II stress intensity component and it will most 
likely experience this kind of oscillation. The explanation of this will be given in the 
following section. 

Relatively small values of the Mode II component of the local stresses make the energy 
release rate due to the crack advance behave similarly to the Mode I stress intensity factor. 
The 7-integral is not path independent in the case of a curvilinear crack (this matter is 
discussed in [2]), and, additionally, due to the presence of the hole, the energy release of the 
system cannot be associated with the local energy release at the crack tip. Thus, the 
traditional thought that the crack will propagate in the direction of maximal energy release 
of the system is not justifiable, and as is shown in Figure 5, the crack may propagate in the 
direction of decreasing energy release rate, as compared with its history along the crack 
path. This means that the decisive factor of the crack growth direction is a strictly local 
preference of the crack tip orientation. 

The data obtained here cannot support the suggestion [11,15] that the derivative 
dKn/dL ( L is the crack length) plays a significant role in the crack curving mechanism. 
However, one may make an argument against the importance of this derivative, since, as the 
data show below, the local variations with respect to a potential deflection angle are more 
important for the crack path formation, than the history type dependence on the crack 
length increment. 

3.2 Variations of the crack tip parameters with respect to possible crack path deflections. 

In this part of the investigation the crack propagation of the numerical model was assumed 
to be partially along the given trajectory up to certain point, and then small perturbations of 
the experimental path were given. These perturbations were given in the form of a circular 
arc tangential to the original path and small in length if compared with any geometrical 
parameter of the problem. The length of these arcs was chosen so small that it did not affect 
the result, and the only important parameter that remained was the angle formed by the 
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tangential at the free end of the arc and by the tangential to the actual crack path trajectory 
at the point correspondent to the increment of the crack advance equal to the arc length. 
These crack path perturbations are illustrated schematically in Figure 6. The length of the 
arcs is enlarged in Figure 6 for better visualization. The sign of the angle of the crack path 
deflection is marked in Figure 6. 

The aim of this analysis is to determine what makes the actual crack path preferable 
with respect to other possible crack path directions. The results are given in Figure la-d for 
three crack tip positions x Q = 2.0, 3.0, 4.0 and deflections of ± 10 degrees. The variation of 
the Mode II stress intensity factor is given for one position only because it exhibits 
practically equivalent behavior for all cases considered. Variations of the Mode I stress 
intensity factor are given for all three cases. The data for the energy release rate are 
omitted here since they practically duplicate the variations of Aij with changed scale. 

The data for K j are normalized by the value of K\(d =0) corresponding to zero 
deflection angle, and the values of Xjj are normalized by the value of the reference Mode I 
stress intensity factor, as in the previous section. 

None of the considered positions can be described as zero K\\ direction, as a direction 
of maximal K\, or as a direction of maximal energy release rate due to the crack advance. By 
the directions of maximal value here, one means maximal value with respect to an angle of 
possible crack path perturbations. Directions of maximal energy release rate due to the 
crack advance are slightly off the course of the crack path. These directions practically 
coincide with K\\ = 0 direction and, therefore, coincide with maximal direction, as it was 
predicted by Cotterell and Rice [9], On Figure lb , directions of zero K\\ are marked on the 
curves. The difference between the maximal ATj directions and the marked directions of 
* 11 = 0 is less than 0.2 degrees, and this is a well reproduced result. The difference between 
the actual crack path direction and direction of a maximal energy release varies between 2 
to 15 degrees. Denoting the deflection angle ce and the energy release rate per crack 
advance G, we observe that for higher value of the derivative dG/da at a =0, the direction 
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of the maximal G is close to the actual path. It seems to be more suitable physically to 
describe this process in terms of the energy release rate rather than in terms of stress 
intensity factor, although quantitatively both would give the same result, since direction of 
maximal K\ coincides with direction of zero /Cjj. 

4. Proposed crack path deflection mechanism. 

To explain the observed data which were described in the previous section, it is 
essential to consider the material properties exhibited during the crack growth process. As 
was pointed out in the introduction, the material which was used in the experiment may be 
characterized by the formation of a narrow nonlinear zone along the crack path. In the 
framework of linear fracture mechanics, the crack and the surrounding zone cannot be 
separated and are treated as a single object. The resulting fracture mechanics parameters 
then may be interpreted as remote, or applied, values controlling the nonlinear region. The 
nonlinear zone acts as a buffer for the crack in absorbing the asymmetry of the applied 
stress field. As a result, the nonlinear zone becomes asymmetric with respect to a crack line. 
The nature of the nonlinear deformation mechanism is not important in principle, and for 
purposes of this study, it may be characterized by the quantity of the stress field asymmetry 
it can absorb before turning the crack. In the framework of linear fracture mechanics, the 
derivative d G/d a at a = 0 and Aqj are the only parameters that may be used as a measure of 
the asymmetry of the applied stress field acting in the vicinity of a crack tip. Of these two 
parameters, the energy release gradient with respect to the deflection angle, dG/da, has a 
better physical interpretation in terms of the energy release rate increase (benefit) due to 
the crack tip rotation. It is a material property that determines the critical value of the 
energy release rate benefit required to rotate the crack tip. On the other hand, considering 
that both parameters play similar roles, one has to note that the critical value of Aiji would 
be easier to measure. In Figure 8, schematic development of an asymmetric crack 
propagation process zone is illustrated. Thus, 


li 



s 


a - b 
a + b 


a 


(3.1) 


8G 

da 


(a=0) 


where b and a are linear dimensions of the wake of the nonlinear process zone, as 
illustrated in Figure 8. The ratio S on the right of (3.1) may be taken as a geometric 
criterion of the asymmetric development of the nonlinear region. For positive values of a 
and b, this ratio changes from -1 to 1, with 0 corresponding to the symmetric case. Sign of S' 
corresponds to the direction of the angle of deflection. 

The critical value of the energy release rate gradient at which material will require a 
change in crack growth direction may determine the instant at which the crack tip 
orientation will change, but it cannot, generally speaking, determine the amount of crack tip 
deflection because of the buffering effect of the nonlinear process zone. The new crack 
growth direction will most likely be dictated by the microstructure of the material at the 
particular instant, such as the grain size, or local orientation of crystallographic planes. At 
this point, the crack will continue to grow in that direction until the critical situation 
develops again, and the energy release rate gradient reaches the critical value for the 
material; consequently, the crack growth direction will change again. Thus, the ductile 
materials with high buffering capability can tolerate a higher degree of stress field 
asymmetry, and require smaller crack path alterations to accommodate it. Therefore, 
ductile materials will exhibit smooth crack path trajectories, and brittle materials with a very 
narrow nonlinear zone or without this buffering zone will exhibit sharp turns along the 
crack growing path. This type of observation was reported in [20], 

In the analyzed experiment, the critical value of the energy release rate gradient (or 
critical value of /^jj) was never reached in case 4, and it was nearly reached in case 3. In 
other cases, as the curvature of the trajectory changes, the crack path experiences a 
sequence of corrections in its direction; this explains the oscillation of K\\ observed earlier. 
The frequency of this oscillation is determined by the material property associated with a 
critical energy release rate gradient and by the intensity of the applied stress field. Case 2 
clearly demonstrates this oscillation of Aqj, as is shown in Figure 4. 
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It follows, that under mixed mode loading, only the ideally brittle materials are capable 
of forming a crack path with Aqj = 0 along the significant crack path segments, and these 
segments will be nearly straight line intervals, 

5. Discussion and conclusions. 

The experimental crack path trajectories were analyzed to establish the governing 
mechanism of the crack path formation. The experimental data were supplied by the 
authors of reference [1], The experimental crack path trajectories were used as boundary 
value data for the numerical analysis. Through accurate numerical procedure the essential 
fracture mechanics parameters were analyzed as they develop along the practical crack path 
trajectories and as they change due to crack path perturbations. 

Results reported here were reproduced from data from different specimens, which 
shows their reliability. 

It was observed that the crack path is not characterized by monotonic increase of the 
energy release rate due to the crack advance, and the crack growth direction is controlled by 
local gradients of the energy release rate with respect to a deflection angle. 

The proposed crack path formation mechanism is based on the assumption that the 
process zone formed during the crack growth serves as a buffer between the applied load 
asymmetry and the crack tip reaction to it. Thus, the magnitude of the energy release rate 
gradient, required to alter the crack growth direction, is associated with the process zone in 
the vicinity of the crack tip. The ability of this zone to transmit the asymmetric properties of 
the applied stress field is a material property, and is important for the crack path alteration. 
In other words, the buffering properties of the process zone must be characterized by the 
material constant. This constant can be associated with the maximal value of the energy 
release rate gradient with respect to deflection angle, which the crack process zone can 
tolerate without a change in crack growth direction. This property of the material may be 
described as a maximal allowable geometrical asymmetry of the process zone, which can be 
specified by the parameter S introduced in equation (3.1). 
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The results of the numerical analysis show that the direction of the maximal energy 
release rate due to the crack advance coincides with the direction of the zero Mode II stress 
intensity factor, which is consistent with [9]. Therefore, the critical value of Aqj, which 
determines the instant when the crack tip changes orientation for a particular material, can 
be used as an equivalent to the energy release rate gradient. The disadvantage of using the 
Kj j criterion is that for the crack growth, the certain value of the Mode I stress intensity 
factor has to be maintained, and the combination of both Modes of the applied stress field 
is important. 

As the crack tip approaches the inhomogeneity, the magnitude of the energy release 
rate gradient increases; after it reaches the critical value for the material, the crack growing 
direction changes until the critical value of the gradient is reached again. If the material is 
homogeneous, this process can be accomplished in small crack growth increments forming a 
smooth crack path curve. As was mentioned, the Mode II stress intensity factor plays a role 
similar to the energy release rate gradient in controlling the crack path. The curvilinear 
crack path experiences continuous correction, and, therefore, the oscillations of K\\ are 
observed along the crack path. Accordingly, only in very brittle materials with negligible 
nonlinear process zone size, the crack path may exhibit sharp turns and may have straight 
line intervals of significant length with near zero values of ATjj. 
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FIGURE 1, - GEOMETRY OF THE SPECIMEN. 


FIGURE 2. - EXPERIMENTAL CRACK PATH TRAJECTORIES (DATA FOR THESE 
TRAJECTORIES WERE SUPPLIED BY AUTHORS 111). 
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KII OF GROWING CRACKS 



FIGURE 3. - Kj VARIATION ALONG THE CURVILINEAR CRACK PATH 
TRAJECTORIES. 



FIGURE - K n VARIATION ALONG THE CURVILINEAR CRACK PATH 
TRAJECTORIES. 



FIGURE 5. - ENERGY RELEASE RATE DUE TO THE CRACK ADVANCE 
ALONG THE CURVILINEAR CRACK PATH TRAJECTORIES. 


TRAJECTORY ANALYSIS 



FIGURE 6. - CRACK PATH DEFLECTIONS IN THE FORM OF CIRCULAR ARCS 
USED IN THE ANALYSIS OF THE LOCAL VARIATIONS OF THE FRACTURE 
MECHANICS PARAMETERS. THE LENGTH OF THE ARCS IS EXAGGERATED 
IN THIS ILLUSTRATION. 
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